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Abstract.  In  this  paper  we  report  our  recent  results  on  the  design  of  materials  with 
controlled  properties  by  the  application  of  computational  chemistry  methods,  new 
algorithms  and  scalable  software. 


1  Introduction 

Developments  of  new  materials  (such  as  new  high  energy  density  materials;  POSS  (polyhedral 
oligomeric  silsesquioxanes)  that  are  highly  resistant  to  extreme  environments  and  are  therefore 
desirable  coatings  for  rocket  engines;  nonlinear  optical  materials;  and  liquid  crystals)  are 
important  in  many  applications  of  interest  to  the  Air  Force.  However,  the  design  of  large  and 
complex  material  systems  with  controlled  properties  is  challenging.  Therefore,  the  calculations 
we  describe  in  this  paper  would  not  have  been  possible  without  scalable  quantum  chemistry 
codes  such  as  GAMESS  (General  Atomic  and  Molecular  Electronic  Structure  System), 
developed  under  the  auspices  of  a  DoD  CHS  SI  (Computational  High  Performance  Scalable 
Software  Initiative)  grant,  Gaussian,  and  molecular  dynamics  codes,  as  well  as  the  availability  of 
the  DoD  Computation  Centers  through  our  computational  challenge  project. 

2  Results  and  Discussion 

2.1  High-Energy  Density  Materials  (HEDM) 


2.1.1  Bond  dissociation  energies:  a  comparative  study  of  theoretical  methodsl’2 

The  purpose  of  this  study  is  to  perform  a  systematic  comparison  of  several  theoretical  methods 
for  computing  the  bond  dissociation  energies  (BDEs)  of  a  prototypical  energetic  compound.  The 

1  This  work  was  done  in  collaboration  with  Donald  L.  Thompson  and  Dan  C.  Sorescu,  Department  of  Chemistry,  Oklahoma  State 
University,  Stillwater,  OK  74078. 

2  We  thank  Dr.  Jeffrey  A.  Sheehy  for  several  stimulating  discussions. 


computational  methods  include  single-reference  (Moller-Plesset  second  order  perturbation  theory 
(MP2)[1],  Gaussian-2  using  MP2  (G2(MP2))[2],  and  coupled  cluster  with  single  and  double 
substitutions  (CCSD(T))[3]),  multiconfigurational  (complete  active  space  self-consistent  field 
(CASSCF)[4]  and  multiconfiguration  quasi-degenerate  perturbation  theory  (MCQDPT)[5]),  as 
well  as  density  functional  theory  (DFT,  using  Becke's  three  parameter  hybrid  exchange 
functional  with  the  Lee,  Yang,  and  Parr  correlation  (B3LYP)[6])  methods.  The  energetic 
compounds  of  interest  are  l,l-diamino-2,2-dinitroethylene  (FOX-7)  and  its  simpler  analogues. 
FOX-7  was  recently  synthesized[7]  and  found  to  have  superior  impact  and  friction  sensitivity 
properties  relative  to  the  well-known  HMX  and  RDX  explosives.  This  molecular  system  is  an 
example  of  a  “push-pull”  ethylene  compound,  which  has  both  electron  donating  (NH2)  and 
electron  withdrawing  (N02)  groups  attached  to  a  common  ethylene  linkage.  Such  compounds  are 
potentially  stabilized  through  “71”  resonance  delocalization  and  hence  may  possess  non- 
negligible  multiconfigurational  character.  Therefore,  this  molecule  should  be  a  good  testbed  for 
comparing  DFT  and  single-configuration  methods  against  multiconfigurational  techniques. 

The  geometries  of  the  reactants  and  products  of  the  C-NH2  and  C-NO2  homolytic  bond 
dissociation  reactions  of  FOX-7  and  its  analogues  have  been  fully  optimized  and  verified  as  local 
minima  at  the  B3LYP,  MP2  and  CASSCF  levels,  using  the  6-311G(d,p)  basis  set. [8] 
CCSD(T)//MP2,  MCQDPT//C AS SCF,  and  G2(MP2)  energies  were  also  computed.  The  putative 
resonance  delocalization  of  the  “71”  electrons  in  FOX-7  suggests  that  in  addition  to  the  (2e,2o)  C- 
N  a, a*  space  needed  to  describe  homolytic  bond  cleavage,  an  appropriate  active  space  for  the 
MCSCF  calculations  should  include  the  (2e,2o)  C-C7t+7t ,  (4e,3o)  N02  7r+non-bonding+7r ,  and 
(2e,lo)  NH2  lone  pair  electrons  and  orbitals,  described  in  further  detail  elsewhere.3  The  degree  of 
multiconfigurational  character  is  assessed  by  calculation  of  the  natural  orbital  occupation 
numbers  (NOONs)  of  the  MP2,  CCSD(T),  and  CASSCF  first  order  density  matrices.  As  shown 
in  a  recent  work  by  Gordon[9]  et  al.,  the  occurrence  of  “non-physical”  MP2  and  CCSD(T) 
NOONs  is  an  indication  that  a  single  configuration  treatment  may  be  inadequate.  All  calculations 
were  performed  using  GAMES S[  10],  GAUSSIAN98[1 1],  and  ACESII.[12] 

At  the  MP2  level,  aminoethylene  is  the  only  closed  shell  compound  without  at  least  one  non¬ 
physical  NOON.  Although  some  calculations  are  still  in  progress,  it  appears  that  all  open-shell 
species  except  the  amino  and  nitro  radicals  have  multiple  non-physical  NOONs  at  this  level  of 
theory.  In  contrast,  none  of  the  compounds  have  non-physical  NOONs  at  the  CCSD(T)  level, 
although  several  open  shell  CCSD(T)  calculations  are  still  in  progress.  This  indicates  that 
CCSD(T)  is  able  to  capture  more  of  the  non-dynamical  correlation  than  MP2  and  hence  should 
be  more  reliable.  The  CASSCF  NOONs  show  that  the  lone  pair  orbital  on  the  NH2  group  has 
little  multiconfigurational  character.  Elowever,  the  nominally  unoccupied  N02  7i*  (0.08-0.14 
electrons)  orbitals  have  relatively  large  populations,  suggesting  that  the  C-N02  BDEs  may  not  be 
as  well  described  by  single  reference  methods  as  the  C-NH2  BDEs. 

The  predicted  C-NH2  BDEs  for  aminoethylene  (which  include  zero-point  energy  corrections)  are 
all  within  9  kcal/mol  of  the  experimental  value  of  102  kcal/mol[13]  and  range  from  93 
(MCQDPT)  to  108  (MP2)  kcal/mol.  The  B3P86/6-31+G(d,p)[14]  and  G2(MP2)  predictions  are 
both  within  about  1  kcal/mol  of  experiment.  MP2  overestimates  the  BDE  by  6  kcal/mol  and 
MCQDPT  underestimates  it  by  about  9  kcal/mol,  with  CCSD(T)  (97.4  kcal/mol)  and  B3LYP 


3  Boatz,  J.A.,  Sorescu,  D.C.,  Thompson,  D.L.:  manuscript  in  preparation 


(96.7  kcal/mol)  yielding  intemediate  values.  CCSD(T)  is  the  single  reference  method  in  closest 
agreement  with  MCQDPT. 

For  l-amino-2-nitroethylene,  MP2  again  gives  the  largest  C-NH2  BDE  (126  kcal/mol)  and 
MCQDPT  the  smallest  (112  kcal/mol).  CCSD(T)  is  again  the  single-reference  method  in  closest 
agreement  with  MCQDPT.  Although  the  C-NH2  BDE  calculations  for  the  entire  set  of  ethylene 
derivatives  are  incomplete,  it  appears  that  (a)  MP2  consistently  gives  the  highest  C-NH2  BDEs 
and  (b)  the  MCQPDT  BDEs  are  in  good  agreement  with  CCSD(T).  In  contrast,  the  DFT  BDEs 
show  no  clear  patterns  of  behavior.  For  nitroethylene,  the  C-N02  BDEs  range  from  59.5 
kcal/mol  (MCQDPT)  to  79.4  kcal/mol  (G2(MP2)).  CCSD(T)  is  the  single-reference  method, 
which  is  in  closest  agreement  with  MCQDPT,  although  the  difference  between  these  two  values 
is  greater  than  in  the  case  of  the  C-NH2  BDEs  of  aminoethylene.  The  B3LYP  value  of  68.8 
kcal/mol  is  in  good  agreement  with  CCSD(T)  (70.9  kcal/mol),  although  the  B3P86  prediction  of 
74.5  kcal/mol[14]  is  closest  to  MP2  (77.7  kcal/mol).  Similar  observations  apply  to  the  predicted 
C-N02  bond  energies  in  l-amino-2-nitroethylene,  which  range  from  69.3  (MCQDPT)  to  88.8 
kcal/mol  (G2(MP2)).  CCSD(T)  (78.2  kcal/mol)  is  the  single-reference  method  in  closest 
agreement  with  the  MCQDPT  result  and  also  is  in  good  agreement  with  the  B3LYP  value  (81.4 
kcal/mol).  MCQDPT  also  predicts  the  smallest  C-N02  BDE  for  l-amino-2-nitroethylene  (69 
kcal/mol).  Although  the  C-N02  BDE  calculations  are  not  yet  complete,  it  appears  that  (a) 
G2(MP2)  and  MP2  consistently  give  the  highest  C-N02  BDEs  and  (b)  the  single  reference 
method  in  best  agreement  with  MCQDPT  is  CCSD(T). 

2.1.2  High-energy  density  materials 

A  recent  major  focus  centers  around  a  series  of  proposed  monopropellants4,  such  as: 


Isp,  which  is  a  measure  of  the  energy  released  by  the  reaction  of  the  fuel  divided  by  the  mass  of 
the  products,  is  of  primary  interest: 


Isp  oc  [AH/m]1/2  (1) 

Using  a  combination  of  isodesmic  reactions  and  the  G2  model,  the  heat  of  formation  for  this 
species  is  predicted  to  be  456.8  kcal/mol.  This  translates  to  an  Isp  of  240  sec,  compared  with  230 
sec  for  hydrazine,  thus  being  a  very  promising  fuel.  However,  it  is  important  to  consider  the 
stability  of  such  high-energy  species  to  various  reactions,  before  asserting  their  viability  as  fuels. 
One  possible  reaction  is  the  loss  of  molecular  nitrogen  to  form  the  smaller  species  shown  below: 


4  Work  of  the  late  Rob  Schmitt  and  Jeff  Botaro,  SRI. 


Using  our  highly  scalable  second  order  perturbation  theory  code  with  a  large  basis  set,  we  have 
found  that  this  decomposition  process  is  exothermic  by  60  kcal/mol.  It  may  be  that  the  barrier  for 
decomposition  is  large,  in  which  case  the  species  may  still  be  viable.  The  search  for  the  transition 
state  for  this  reaction  is  in  progress. 

Recently  the  first  new  all-nitrogen  compound  in  nearly  a  century,  Ns+[15],  was  synthesized.  The 
counterion  for  this  species  was  AsF6’,  the  result  being  a  rather  unstable  solid.  In  order  to 
understand  the  potential  energy  surface  for  this  species,  we  performed  MP2  calculations  on  the 
various  stationary  points.  The  first  step  in  this  process  was  to  separate  the  two  ions  by  a  large 
distance  and  optimize  the  geometry.  The  result  was  an  ion  pair  [AsF6'][Ns+]  complex  that  is,  not 
surprisingly,  98  kcal/mol  below  the  separated  ions.  What  was  somewhat  surprising  is  that  this 
complex  is  47  kcal/mol  below  the  neutral  species  that  is  formed  when  one  F  ion  is  transferred  to 
N5+  to  form  neutral  AsFs  +  FN5.  The  former  is  a  well-known  species,  but  the  latter  has  not 
previously  been  reported.  FN5  has  several  isomers,  whose  relative  energies  and  the  barriers 
separating  them  were  determined  using  large  basis  set  CCSD(T)  calculations.  The  lowest  lying 
isomers  were  found  to  have  small  energy  barriers  (<  10  kcal/mol)  separating  them  from  the  much 
more  stable  decomposition  products  FN3  +  N2.  The  exothermicity  of  this  decomposition  is  more 
than  50  kcal/mol,  so  the  products  AsFs  +  FN3  +  N2  are  slightly  lower  in  energy  than  the  original 
ion  pair  [AsF6"][Ns+].  Since  FN3  is  itself  an  energetic  species  with  a  relatively  small  barrier  to 
further  decomposition,  it  is  therefore  not  surprising  that  the  original  species  made  by  Christe  et 
al.  is  unstable. 

2.1.3  Quantum  molecular  dynamics  (MD)  simulations  of  low  temperature  HEDM 

Traditional  MD  algorithms  can  be  extended  to  incorporate  equilibrium  quantum  mechanical 
effects  through  the  use  of  discretized  Feynman  path  integrals [16],  The  resulting  algorithm,  path 
integral  MD  (PIMD),  is  highly  amenable  to  parallelization  because  each  physical  particle 
becomes  rigorously  mapped  onto  a  collection  of  classical-like  quasiparticles  at  different  values 
of  imaginary  time,  each  connected  to  its  neighbor  by  a  harmonic-like  force.  The  discretized 
PIMD  algorithm  can  be  load-balanced  across  computational  nodes  because  the  quasiparticles 
with  the  same  imaginary  time  index  experience  identical  forces.  The  inter-node  communications 
are  minimal  since  there  are  only  nearest-neighbor  harmonic  forces  to  compute.  The  PIMD 
algorithm,  no  matter  how  complex  the  physical  system,  can  always  be  made  to  achieve  excellent 
scalability  and  thus  high  performance.  Moreover,  classical  MD  techniques  such  as  constant 
temperature  and  pressure  algorithms  can  be  readily  incorporated  into  PIMD.  The  algorithm  is 
even  amenable  to  a  two-tier  level  of  parallelism,  the  first  being  over  the  imaginary  time  slices 
(quasiparticles),  the  second  over  the  force  loop  for  the  interparticle  interactions  in  large 
(thousand  or  even  million)  particle  simulations. 

In  an  important  theoretical  advance[17,18,19,20],  PIMD  methods  have  been  extended  to 
include  quantum  dynamical  effects  (i.e.,  not  just  equilibrium  properties).  This  computational 


approach  is  called  "Centroid  Molecular  Dynamics"  (CMD).  It  incorporates  the  dominant 
quantum  dynamical  effects  of  a  many-body  system  into  a  classical-like  MD  framework,  thus 
significantly  extending  the  range  and  applicability  of  MD  methods.  The  basic  result  of  CMD  [19- 
22]  is  that  the  dynamical  correlations  of  quantum  particles  in  a  general  many-body  system  can  be 
accurately  computed  by  running  classical-like  trajectories  on  an  effective  potential  which 
includes  the  effects  of  quantum  zero  point  energy  and  tunneling.  The  task  of  integrating  the 
CMD  equations  is  not  trivial  since  the  effective  potential  is  a  quantum  potential  of  mean  force, 
requiring  an  "on  the  fly"  dynamical  averaging  at  each  timestep.  A  number  of  powerful 
algorithms  for  solving  the  CMD  equation  have  been  developed[22].  The  best  CMD  algorithm  so 
far  is  the  “hyperparallel”  CMD  approach[21]  using  an  adiabatic  algorithm[22]  in  the  first  tier  of 
parallelism  to  determine  the  effective  potential;  a  second  tier  of  parallelism  can  be  introduced  in 
which,  for  each  of  the  averaging  slaves,  the  MD  steps  are  sent  to  additional  nodes  for  each 
discretized  path  integral  imaginary  time  slice[23],  A  third  tier  of  parallelism  for  very  large 
systems  can  be  utilized  for  each  imaginary  time  slice,  in  which  the  classical  force  loop  is  split 
over  several  nodes.  The  scaling  of  the  hyperparallel  CMD  code  on  the  IBM  SP  and  SGI  02000, 
as  well  as  a  Pentium  Linux  cluster,  is  excellent  (Fig.  1).  The  parallelization  of  the  long  rage 
forces  (Ewald  summation)  was  also  achieved. 


Fig.  1.  Scaling  of  CMD 

The  PIMD  and  CMD  algorithms  are  of  general  interest  in  all  areas  of  condensed  matter 
computer  simulation.  However,  they  are  particularly  important  for  the  HEDM  program  in  the 
computational  study  of,  e.g.,  solid  hydrogen,  helium,  or  nitrogen  matrices  which  contain  isolated 
reactive  species  such  as  lithium,  boron,  aluminum,  hydrogen  atoms,  organic  radicals,  etc.  Such 
matrices  are  high  priority  targets  as  possible  HEDMs  for  Air  Force  space  propulsion  purposes. 
At  the  temperatures  appropriate  to  the  condensed  phases  of  hydrogen  or  helium,  the  host  matrix 
molecules  will  exhibit  large  quantum  mechanical  effects,  but  the  treatment  of  these  effects  by  a 
direct  attack  on  the  time-dependent  Schrodinger  equation  is  impossible.  The  key  issue  is  the 
stability  of  the  HEDM  (i.e.,  the  ability  to  trap  the  energetic  impurities  for  some  period  of  time 
and  thus  impede  their  diffusion  and  eventual  recombination).  Both  the  structural  and  dynamical 
aspects  of  the  composite  condensed  phase  HEDM  systems  are  of  interest.  For  the  equilibrium 
structural  studies,  PIMD  methods  are  adequate.  For  the  dynamical  studies,  the  CMD  approach 
must  be  employed  to  directly  simulate  the  diffusive  and  recombination  steps  of  the  guest  atoms 
in  the  solid  host.  The  long-term  goal  of  this  research  is  to  develop  general  computational 
methods  to  rapidly  and  efficiently  characterize  proposed  high  energy  density  materials. 

The  structure  and  stability  of  atomic  impurities  trapped  in  solid  />ara-hydrogen  have  been 
studied  by  employing  large  scale  PIMD  and  CMD  simulations  at  4  K  and  zero  external  pressure. 
Starting  from  pure  solid  hydrogen,  doped  systems  were  prepared  by  substituting  impurity  atoms 
for  hydrogen  molecules  at  substitutional  defect  sites.  The  structural  and  thermodynamic 


properties  of  the  resulting  HEDM  were  then  calculated.  In  the  case  of  atomic  boron  impurities, 
the  impurity  atom  interacts  anisotropically  with  the  matrix  _para-hydrogen  molecules  because  of 
its  singly  filled  2p  orbital.  To  assess  the  effect  of  the  electronic  anisotropy  in  the  potential,  a 
comparison  has  been  made  between  simulations  in  which  the  orientation-dependent  (anisotropic) 
and  orientation-averaged  (spherical)  B— H2  potentials  were  used.  In  collaboration  with  the 
Alexander  group[22],  three  matrices  were  investigated:  (a)  one  with  a  single  B  atom  site- 
substituted  for  a  p-Ho  molecule,  (b)  one  with  a  similar  site-substituted  matrix  having  a  nearest- 
neighbor  vacancy,  (c)  and  one  with  a  periodic  "super  cell"  simulation  box  having  two  hydrogen 
bulk  slab  regions  in  order  to  study  a  B  impurity  near  the  p-H.2  surface. 

It  was  found  that  significant  “second-order”  distortions  of  the  lattice  occur  to  permit  an 
energetically  favorable  orientation  of  the  2p  orbital,  even  in  the  absence  of  a  neighboring 
vacancy.  Especially  when  the  B  impurity  is  located  near  the  surface,  the  spherically  averaged 
potential  provides  a  noticeably  different  description  from  the  case  of  the  anisotropic  potential. 
These  quantum  MD  simulations  reveal  the  large  structural  differences  that  can  be  found  in  the 
lattice  resulting  from  the  anisotropic  boron-hydrogen  interaction  potential.  The  implications  of 
this  behavior  will  be  explored  in  future  research. 

2.2  POSS  Compounds 

There  is  great  interest  in  POSS  compounds  because  of  their  resistance  to  extreme  environments 
and  the  ease  with  which  they  are  synthesized.  Since  almost  nothing  is  known  experimentally 
about  the  mechanism  by  which  they  are  formed,  we  have  embarked  on  a  long  term  project  to 
determine  the  possible  competing  mechanisms  as  a  function  of  catalyst,  solvent,  and  substituents 
(e.g.,  in  [23]).  The  synthesis  begins  with  hydrolysis  of  RS1X3  to  RSi(OH)3,  followed  by 
condensation  to  the  siloxane  RSi(0H)2-0-Si(0H)2R.  Subsequent  condensations  lead  to  the  3D 
cage  compounds.  We  have  shown  that  the  initial  hydrolysis  and  condensation  steps  all  have  very 
high-energy  barriers  in  the  gas  phase,  but  are  being  reduced  to  nearly  zero  by  the  presence  of  one 
water  molecule  added  to  represent  the  solvent.  Thus,  all  steps  leading  to  the  initial  disiloxanes 
and  to  the  ring  compounds  D3  and  D4  occur  with  net  energy  requirements  of  less  than  10 
kcal/mol.  In  this  paper  we  describe  current  efforts  to  (a)  determine  substituent  effects  on  these 
barriers,  (b)  determine  the  effects  of  adding  additional  water  molecules,  and  (c)  to  compare  the 
properties  of  the  3D  POSS  compounds  with  their  Ti  analogs.  The  substituents  studied  include 
R=H,  methyl,  t-butyl  and  phenyl,  and  X=C1,  OCH3.  There  is  experimental  interest  in 
synthesizing  incompletely  condensed  POSS.  Because  the  reaction  is  very  fast,  it  always  goes  to 
completion.  The  interest  in  substituent  effects  lies  in  attempting  to  slow  down  the  reaction, 
possibly  using  bulky  R  groups.  There  is  industrial  interest  in  the  possibility  that  Ti  POSS  may  be 
another  class  of  materials  with  desirable  properties.  Hammes-Schiffer  et  al.  [24]  have  been 
developing  methodologies  for  the  simulation  of  hydrogen  transfer  reactions,  and  most  recently  a 
new  molecular  orbital  method  for  the  simultaneous  calculation  of  electronic  and  hydrogen 
vibrational  wavefimctions,  now  incorporated  into  GAMES  S.  This  method  may  be  used  to  obtain 
minimum  energy  reaction  paths  and  direct  dynamics  trajectories  of  hydrogen  transfer  reactions, 
with  the  advantage  that  nuclear  quantum  effects  such  as  zero  point  energy  and  hydrogen 
tunneling  are  incorporated  during  the  generation  of  the  reaction  paths  and  trajectories,  rather  than 
subsequently  included  as  corrections.  This  approach  will  be  combined  with  a  mixed 
quantum/classical  surface  hopping  method  to  study  the  quantum  dynamics  of  hydrogen  transfer 
reactions  with  an  ah  initio  potential  energy  surface  obtained  “on  the  fly,”  and  used  to  investigate 


hydrogen  transfer  reactions  in  the  hydrolysis  and  condensation  steps  required  for  the  synthesis  of 
POSS,  in  order  to  determine  the  influence  of  different  trihalosilane  reactants  on  the  rates  and 
yields  of  these  steps  and  aid  in  the  efficient  synthesis  of  POSS. 


2.3  Nonlinear  Absorbing  Materials 

We  have  previously  shown  that  density  functional  theory  and  the  time-dependent  density 
functional  response  theory  (TDDFT)  (B3LYP)  accurately  predict  the  structures  and  spectra  of 
porphin  and  the  larger  tetraphenylporphyrin  and  its  (3-octahalogenated  derivatives, [25,26]  thus 
enabling  their  assessment  as  reverse  saturable  absorbing  materials.  For  porphyn  (PFQ,  the 
photoinduced  tautomerization  process  has  also  been  investigated[27]  in  order  to  discern  between 
previously  proposed  discordant  mechanisms.  Recently  we  carried  out  calculations[28]  for 
tetrabenzoporphyrin  (ZnTBP)  and  phthalocyanine  (ZnPc)  (Fig.  2). 


Fig.  2.  Zinc  phthalocyanine 

Our  calculations  support  the  multiple  bands  in  the  B  and  higher  energy  regions  of  the  ZnPc 
and  ZnTBP  spectra  that  have  been  determined  from  experiment.  Despite  considerable 
experimental  and  theoretical  work  on  ZnPc  and  ZnTBP,  their  ground  state  electronic  structures 
and  UV  spectra  are  not  well  understood.  Previous  theoretical  studies  provided  only  a  qualitative 
understanding  in  the  low-energy  region,  and  none  for  the  higher  energy  bands.  We  have 
elucidated  significant  shifts  in  Gouterman’s  four  orbitals  of  ZnP  upon  tetraazasubstitution  and 
tetrabenzoannulation.5  The  benzo  groups  destabilize  the  lalu  orbital  in  ZnP,  but  they  do  not  show 
a  significant  deviation  from  the  four-orbital  model  for  the  interpretation  of  the  Q  and  B  bands  of 
ZnTBP.  However,  tetraazasubstitutions  considerably  lower  the  4a2u  orbital  (HOMO-1)  in 
ZnTBP.  Also,  it  is  shown  that  the  near  degeneracy  of  the  HOMO  and  HOMO-1  in  ZnP  that 
provides  the  basis  for  Gouterman’s  four  orbital  model  does  not  hold  for  ZnPc.  In  ZnPc,  the  near 
degeneracy  of  the  HOMO-1  (4a2u)  with  other  occupied  orbitals  gives  rise  to  a  complex  structure 
of  the  higher  energy  regions  of  the  spectrum  that  were  not  reproduced  by  simple  MO  methods. 
The  4a2u  MO  has  significant  k  contributions  from  the  benzo  rings  in  ZnPc  but  not  in  ZnTBP.  The 
Q  and  B  bands  are  predicted  to  appear  at  2.18  and  3.28  eV  for  ZnTBP,  respectively,  in  excellent 
agreement  with  the  maxima  located  at  2.02  and  3.18  eV  as  established  from  the  absorption 
spectrum  in  a  supersonic  jet  expansion,  with  similar  values  (2.02  &  3.10  eV)  obtained  from  the 
absorption  spectrum  in  an  Ar  matrix.  Tetraazasubstitutions  are  predicted  to  redshift  the  Q  band 
of  ZnTBP  by  about  0. 1  eV.  The  computed  excitation  energies  in  near  quantitative  agreement 
with  experiment  clearly  provide  a  strong  basis  for  the  interpretation  of  the  electronic  spectra  of 
such  large  molecular  systems.  TDDFT  calculations  are  also  being  carried  out  for  candidate  two- 
photon  absorption  (TP A)  materials,6  to  be  compared  to  those  of  Albota  et  al.  [29] 


5  Nguyen,  K.A.  and  Pachter,  R.  J.  Chem.  Phys.,  manuscript  in  preparation 

6  Day,  P.N.  and  Pachter,  R.,  work  in  progress 


2.4  Liquid  Crystals 

Recent  work  is  in  progress  for  studying  a  liquid  crystalline  droplet  to  model  the  behavior  of  a 
liquid  crystal  in  the  bulk,  with  the  details  reported  elsewhere.  [30]  Indeed,  atomic  level 
simulations  are  proven  important  in  understanding  the  structure-property  relations  of  materials. 
Due  to  the  size  of  the  LC  droplet,  very  large-scale  classical  MD  simulations  of  such  systems  are 
required,  that  may,  become  computationally  prohibitive.  The  fast  multipole  method  (FMM)[31], 
which  uses  a  multiscale  hierarchy  of  partitions  of  the  volume  and  a  divide-and-conquer  strategy 
to  compute  the  power  series,  allows  all  the  forces  to  be  computed  to  any  specified  accuracy  in 
<9(N)  operations.  FMM3D  (an  implementation  and  improvement  to  FMM  in  3-D),  which 
contains  a  variety  of  schemes  for  computing  multipole  translations,  has  been  implemented  into 
the  MD  program  that  is  being  developed  in  our  group. [32] 
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of  formation,  high  densities 

-  role  of  hydrogen  bonding  in  reduced  sensitivity,  importance  of  non-additive 
inter/intramolecular  interactions,  decomposition  mechanisms 


xanes  --  next  generation  plastics 
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IP  =  Ionization  Potential,  S0  -Sn=  Ground  State  Spectrum, 

S0  -Tj  =  Singlet-Triplet  Gap  Tt  -Tn  =  Triplet-Triplet  Spectrum 

E  =  Experiment,  C  =  Calculated 
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Development  of  a  solid-state  force  field:  calibration  of  theoretical  methods  for  computing  bond 
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